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The interstellar medium of the Milky Way is multi- 
phase,^ magnetized^ and turbulent.^ Turbulence 
in the interstellar medium produces a global cas- 
cade of random gas motions, spanning scales rang- 
ing from 100 parsecs to 1000 kilometres.^ Funda- 
mental parameters of interstellar turbulence such as 
the sonic Mach number (the speed of sound) have 
been difficult to determine because observations 
have lacked the sensitivity and resolution to directly 
image the small-scale structure associated with tur- 
bulent motion.^ ^ Observations of linear polariza- 
tion and Faraday rotation in radio emission from the 
Milky Way have identified unusual polarized struc- 
tures that often have no counterparts in the total 
radiation intensity or at other wavelengths,^^^^ and 
whose physical significance has been unclear .^^^^^ 
Here we report that the gradient of the Stokes vec- 
tor (Q, LT), where Q and C/are parameters describing 
the polarization state of radiation, provides an im- 
age of magnetized turbulence in diff'use ionized gas, 
manifested as a complex filamentary web of discon- 
tinuities in gas density and magnetic field. Through 
comparison with simulations, we demonstrate that 
turbulence in the warm ionized medium has a rela- 
tively low sonic Mach number, Ms ^ 2. The develop- 
ment of statistical tools for the analysis of polariza- 
tion gradients will allow accurate determinations of 
the Mach number, Reynolds number and magnetic 
field strength in interstellar turbulence over a wide 
range of conditions. 

We consider radio continuum images of an 18-square- 
degree patch of the Galactic plane, observed with the 
Australia Telescope Compact Array (ATCA) at a frequency 
of 1.4 GHz. Data were simultaneously recorded in total 
intensity (Stokes parameter /) and in linear polarization 
(Stokes parameters Q and U). The Stokes / image (Fig. 1) 
shows a typical distribution of radio emission, consisting of 
supernova- remnant shells, ionized regions around massive 
stars (Hii regions), and unresolved distant radio sources. 
However, the corresponding images of Q, [/, and the lin- 
early polarized intensity P = (Q^ + [/^)^/^ in Figure 1 are 
filled with complex structure that bears little resemblance 
to the Stokes / image, as has also been seen in many other 



polarimetric observations at radio frequencies .'^'^^ The in- 
tensity variations seen in Q, [/ and P are the result of small- 
scale angular structure in the Faraday rotation induced by 
ionized gas,^ and are thus an indirect representation of tur- 
bulent fluctuations in the free-electron density and magnetic 
fleld throughout the interstellar medium (ISM). 

A limitation of previous studies is that they usually in- 
terpreted the data in terms of the amplitude, P, and/or 
angle, = | tan~^([//Q), of the complex Stokes vector 
P = (Q^U). However, neither polarization amplitude nor 
polarization angle is preserved under arbitrary translations 
and rotations in the Q — U plane. These can result from one 
or more of a smooth distribution of intervening polarized 
emission, a uniform screen of foreground Faraday rotation, 
and the effects of missing large-scale structure in an inter- 
ferometric data set. In the most general case, we are thus 
forced to conclude that the observed values of P and do not 
have any physical signiflcance, and that only measurements 
of quantities that are both translationally and rotationally 
invariant in the Q — U plane can provide insight into the 
physical conditions that produce the observed polarization 
distribution. 

The simplest such quantity is the spatial gradient of P, 
i.e., the rate at which the polarization vector traces out a 
trajectory in the Q — U plane as a function of position on the 
sky. The magnitude of the gradient is unaffected by rotation 
and translation, and so has the potential to reveal proper- 
ties of the polarization distribution that might otherwise be 
hidden by excess foreground emission or Faraday rotation, 
or in data sets from which large-scale structure is missing 
(as is the case for the data shown in Fig. 1). The magnitude 
of the polarization gradient is: 




The expression in Equation ([T]) can be calculated simply, and 
the corresponding image of |VP| (Fig. 2) reveals a complex 
network of tangled fllaments. In particular, all regions in 
which |VP| is high consist of elongated, narrow structures 
rather than extended patches. In the inset to Figure 2, we 
plot the direction of VP for a small subregion of the image. 
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Figure 1. | Total intensity (i) and linearly polarized intensity (Q, U, P) for an 18-deg^ region of the Southern Galactic 
Plane Survey.l^^ All four images have been generated ■'^■'^ from a set of observations -"^^ taken at the ATCA over the period 1997 April 
to 1998 April using a 96-MHz bandwidth centred on an observing frequency of 1384 MHz. The field is a mosaic of 190 pointings each 
with a total integration time of 20 min, resulting in an approximately uniform sensitivity over most of the field, of 0.8 mJy per beam 
(Stokes /) or 0.55 mJy per beam (Stokes Q & ?7) at an angular resolution of 75 arcseconds (1 Jy = 10~^^ W m~^ Hz"-*^). The scale for 
each image is shown on the right of each panel. The Stokes / image is displayed over a range of -40 to +150 mJy per beam. Because 
the ATCA is an interferometer, it is not sensitive to structure on angular scales larger than 35 arcminutes. Faint wisps can be seen, 
corresponding to the sharp edges of large-scale structures. However, the bulk of the smooth radio emission from Galactic cosmic rays is 
not detected. Imaging artefacts in the form of grating rings and radial streaks can be seen around a few very bright sources, but these 
regions were not used in our statistical analysis. The Stokes Q and U images are displayed over a range of -15 to +15 mJy per beam, 
and the P image covers a range of to 15 mJy per beam. Almost none of the structure seen in Q, U and P has any correspondence with 
any emission seen in Stokes /; the mottled structure results from spatial fiuct nations in Faraday rotation in the ISM. 



demonstrating that VP changes most rapidly along direc- 
tions oriented perpendicular to the elongation of the fila- 
ments. We can explore the frequency dependence of these 
filaments, because the 1.4-GHz ATCA data shown in Fig- 
ure 1 consist of nine independent spectral channels of width 
8 MHz, spread over a total bandwidth of 96 MHz. We 
have constructed images of |VP| for each individual spec- 
tral channel, and these show the same set of specific features 
as in the overall image, albeit at reduced signal-to-noise ra- 
tios. The lack of frequency dependence indicates that the 
high-gradient structures seen in this data set correspond to 
physical features in the ISM rather than to contour lines 
introduced by the particular combination of observing fre- 
quency and angular resolution used.^^^^ 

We first consider the possibility that these filaments of 
high gradient are intrinsic to the source of emission. Abrupt 
spatial transitions in the strength or geometry of the mag- 
netic field in a synchrotron-emitting region would gener- 
ate a large gradient in {Q,U). However, processes of that 
sort would also produce structure in the overall synchrotron 
emissivity, such that we would observe features in the image 
of Stokes / that match those seen in |VP|. No such cor- 
respondence is observed, demonstrating that the regions of 
high polarization gradient are not intrinsic to the source of 



polarized emission but must be induced by Faraday rotation 
in magneto-ionized gas. 

Because the amount of Faraday rotation is proportional 
to the line integral of neB\\ from the source to the observer 
(where Ue is the density of free electrons and B\\ is the uni- 
form component of the line-of-sight magnetic field) , the fila- 
mentary structure seen in |VP| must correspond to bound- 
aries across which Ue and/or B\\ show a sudden increase 
or decrease over a small spatial interval. Such discontinu- 
ities could be shock fronts or ionization fronts from dis- 
crete sources, as have been observed in polarization around 
the rims of supernova remnants, Hii regions and planetary 
nebulae. We have examined this possibility by carefully 
comparing our image of |VP| with images and gradient im- 
ages of Stokes / (tracing shock waves seen in synchrotron 
emission),^ 2 1-cm Hi emission (tracing atomic hydrogen)'^ 
and 656.3-nm Ha (tracing ionized hydrogen) over the 
same field, but do not find any correspondences. 

We conclude that the features seen in |VP| are a 
generic component of diffuse, ionized gas in this direction 
in the sky. To test this hypothesis, we performed a series of 
three-dimensional isothermal simulations of magnetohydro- 
dynamic turbulence in the ISM, each with different param- 
eters for the sonic Mach number, defined as A4s = (|v|/cs), 
where v is the local velocity, Cs is the sound speed, and 
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Figure 2. | The gradient image of linear polarization, |VP|, for an 18-deg^ region of the Southern Galactic Plane 
Survey. |VP| has been derived by applying Equation ([T]) to the Q and U images from Fig. 1; note that |VP| cannot be constructed 
from the scalar quantity P = (Q^ + t/^)-*^/^, but is derived from the vector field P = (Q,U). |VP| is a gradient in one dimen- 
sion, for which the appropriate units are (beam)~^-^. Because P measures linearly polarized intensity in units of millijanskys per 
beam, |VP| has units of millijanskys per (beam)-*^-^. The scale showing |VP| is shown to the right of the image, and ranges from 
to 15 mJy per (beam)-*^-^. The inset shows an expanded version of the structure with highest |VP|, covering a box of side 0?9 cen- 
tred on Galactic longitude 329? 8 and Galactic latitude +1?0. Plotted in the inset is the direction of VP at each position, defined as 
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the amplitude of the gradient is greater than 5 mJy per (beam)-*^'^. 



For clarity, vectors are only shown at points where 



the averaging (indicated by angle brackets) is done over 
the whole simulation. For each simulation, we propagated 
a uniform source of polarized emission through the distri- 
bution of turbulent, magnetized gas. The resultant Faraday 
rotation produces a complicated distribution on the sky of 
Stokes Q and [/, from which we generated a map of the po- 
larization gradient using Equation ([T]). Images of |VP| for 
representative simulations of the subsonic, transonic and su- 
personic regimes are shown in Figure 3. Narrow, elongated 
filaments of high polarization gradient are apparent in each 
simulation in Figure 3, although they differ in their mor- 
phology and degree of organization. In particular, the su- 
personic case (Fig. 3c) shows localized groupings of very 
high-gradient filaments, corresponding to ensembles of in- 
tersecting shocks.^ '^^ By contrast, the subsonic (Fig. 3a) 
and transonic (Fig. 3b) cases show more-diffuse networks of 
filaments, representing the cusps and discont inuities charac- 
teristic of any turbulent velocity field .^^^^^^^ 

Visual comparison of the simulated distributions of 
|VP| with real data (Fig. 2) suggests that the subsonic and 
transonic cases shown in Figures 3a, b more closely resem- 
ble the observations than does the supersonic case. We can 
quantify this statement by calculating the third-order mo- 
ment (skew, 7) and the fourth-order (kurtosis, (3) moments 
of the probability distribution function of |VP| for both 
observations and simulations: these quantities parameterize 
the degree of Gaussian asymmetry in the probability distri- 



bution function, and hence provide information on the level 
of compression due to shocks in the data. ^ 

In the simulations, we found that both the skew and 
kurtosis of |VP| increased monotonically with sonic Mach 
number. We used a genetic algorithm to determine that 
the threshold for strongly supersonic turbulence was 7 > 
1.1 and (3 > 1.5. We then computed the third- and fourth- 
order moments for the observed distribution of |VP| shown 
in Figure 2, and found that 7 = 0.3 and /S = 0.9. 

This analysis of the moments of the polarization gra- 
dient therefore confirms quantitatively what we concluded 
above from visual inspection: the turbulent, ionized ISM in 
this direction is subsonic or transonic. The findings we ob- 
tained by imaging the polarization gradients produced by in- 
terstellar turbulence are supported by recent statistical stud- 
ies of Ha emission measures and of 21-cm Hi column densi- 
ties over large volumes, which hav e sim ilarly found A^s ^ 2 
for warm gas throughout the ISM.'^^^ 

In the simulations shown in Figure 3, the sharp gradi- 
ents in {Q,U) occur as a result of localized high values of 
the gas density and magnetic field, resulting from vortic- 
ity or shock compression. However, the filamentary features 
seen in |VP| may not be easily observable in other types of 
data: for example, if we adopt typical parameters for warm 
ionized gas^^^ of rie ~ 0.3 cm~^ and B\\ ^ 2 /xG, even 
the compression associated with a strong adiabatic shock 
produces across-filament changes in emission measure and 
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Figure 3. The gradient image, |VP|, derived from propagation of linear radio polarization through three different 
isothermal simulations of magnetized turbulence. Each simulation is a 512 x 512 x 512-element periodic box with a linear 
dimension of 0.15 parsecs for each pixel, evolved in time using an essentially nonoscillatory scheme .^^^^^ Three such simulations are 
shown, each labelled with its corresponding value for A4s — subsonic (A4s < 1; a), transonic {A4s ~ 1; b) and supersonic cases 
{A4s > 1; c). At the start of each simulation, the electron density had a uniform value Ue = 0.1 cm~^ and the magnetic field was 
oriented in the plane of the sky with a uniform amplitude of B = 0.3 /iG (subsonic), 1 fiG (transonic) or 2 /iG (supersonic), which 
corresponds to a constant Alfvenic Mach number of A4a = 2 in each case. Turbulence was driven solenoidally in Fourier space at large 
scales (small wavenumber) until the turbulent cascade had fully developed and a steady state between input energy and dissipation had 
been reached. In each case, we illuminated the simulation volume with a background radio source of uniform polarization at an emission 
frequency of 1.4 GHz, with Q = 100 mJy per pixel and U = mJy per pixel at all positions. At each pixel, the line integral of ne^y was 
computed, and the corresponding Faraday rotation applied to the polarized signal, to calculate values of Q and U. No effects of finite 
angular resolution, depolarization or incomplete interferometric visibility coverage were included, so that the observed polarized signal is 
P = 100 mJy per pixel at all positions. We then calculated the gradient, |VP|, using Equation (^. The scales showing |VP| are shown 
on the right of the images, and range from to 25 mJy per (pixel) -"^-^ for panel a, to 100 mJy per (pixel) -"^-^ for panel b, and to 500 
mJy per (pixel) -"^'^ for panel c. 



Faraday rotation measure of only 0.5 parsecs cm~ and 
^ 5 radians m~^, respectively, assuming a spatial scale ^'-^ 
for these structures of ~ 0.5 parsecs. This is below observ- 
able levels in Ha and other tracers of emission measure. The 
rotation measure gradient^ across these interfaces is poten- 
tially observable in spectropolarimetric radio data, but the 
addition of single-dish observations is required to recover 
the total power of the polarized signal. By contrast, even 
a small gradient in rotation measure can produce an arbi- 
trarily large value of |VP| (irrespective of whether single- 
dish measurements are present in the data), provided that 
there is a strong source of background polarized emission 
through which the discontinuities in Faraday rotation are 
viewed. Further investigation of the polarization gradient 
and its statistical properties will provide robust estimates of 
poorly constrained parameters of turbulent flows such as the 
sonic and Alfvenic Mach numbers, the characteristic mag- 
netic field strength, the Reynolds number and the physical 
scale of energy injection. 
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